library(slingshot)
## Loading required package: princurve
library(RColorBrewer)
library(mgcv)
## Loading required package: nlme
## This is mgcv 1.8-23. For overview type 'help("mgcv-package")'.
library(tradeR)
library(edgeR)
## Warning: package 'edgeR' was built under R version 3.5.1
## Loading required package: limma
## Warning: package 'limma' was built under R version 3.5.1
## 
## Attaching package: 'limma'
## The following object is masked from 'package:BiocGenerics':
## 
##     plotMA
library(rafalib)
library(wesanderson)

data <- readRDS("~/Dropbox/PhD/Research/singleCell/trajectoryInference/trajectoryDE/tradeRPaper/simulation/sim2_dyntoy_bifurcating_4/bifurcating_4.rds")
counts <- t(data$counts)
falseGenes <- data$feature_info$feature_id[!data$feature_info$housekeeping]
# this dataset has no null genes.
set.seed(5)
null1 <- t(apply(counts,1,sample))
dimnames(null1) <- list(paste0("H",1:501),paste0("C",1:2011))
null2 <- t(apply(counts,1,sample))
dimnames(null2) <- list(paste0("H",502:1002),paste0("C",1:2011))
null3 <- t(apply(counts,1,sample))
dimnames(null3) <- list(paste0("H",1003:1503),paste0("C",1:2011))
counts <- rbind(counts,null1,null2,null3)
nullGenes <- rownames(counts)[substr(rownames(counts),1,1)=="H"]

pal <- wes_palette("Zissou1", 12, type = "continuous")
truePseudotime <- data$prior_information$timecourse_continuous[colnames(counts)]
g <- Hmisc::cut2(truePseudotime,g=12)

# quantile normalization
FQnorm <- function(counts){
  rk <- apply(counts,2,rank,ties.method='min')
  counts.sort <- apply(counts,2,sort)
  refdist <- apply(counts.sort,1,median)
  norm <- apply(rk,2,function(r){ refdist[r] })
  rownames(norm) <- rownames(counts)
  return(norm)
}
normCounts <- FQnorm(counts)

## dim red
pca <- prcomp(log1p(t(normCounts)), scale. = FALSE)
rd <- pca$x[,1:3]
plot(rd, pch=16, asp = 1)

set.seed(9)
cl <- kmeans(rd, centers = 8)$cluster
plot(rd, col = brewer.pal(9,"Set1")[cl], pch=16, asp = 1)
legend("topleft",legend=as.character(1:7),col=brewer.pal(9,"Set1")[1:7],pch=16,cex=2/3,bty='n')

 plot(rd, col = pal[g], pch=16, asp = 1)
#lineages
lin <- getLineages(rd, cl, start.clus=6, end.clus=c(5,2))
## Using full covariance matrix
plot(rd, col = pal[g], pch=16, asp = 1)
lines(lin,lwd=2)

#curves
crv <- getCurves(lin)
plot(rd, col = pal[g], pch=16, asp = 1)
lines(crv, lwd=2, col="black")

# milestone ID
gid <- data$prior_information$groups_id
gid <- gid[match(colnames(counts),gid$cell_id),]
plot(rd, col = as.numeric(as.factor(gid$group_id))+1, pch=16, asp = 1)

fit smoothers on raw data

cWeights <- slingCurveWeights(crv)
pseudoT <- slingPseudotime(crv, na=FALSE)
gamList <- fitGAM(counts, pseudotime=pseudoT, cellWeights=cWeights, verbose=FALSE)

Test for differences at end point

endRes <- diffEndTest(gamList)
hist(endRes$pvalue)

deGenesEndGam <- rownames(counts)[which(p.adjust(endRes$pvalue,"fdr")<=0.05)]
mean(falseGenes%in%deGenesEndGam) #TPR
## [1] 0.6626747
mean(!deGenesEndGam%in%falseGenes) #FDR
## [1] 0

edgeR analysis on final clusters

clF <- as.factor(cl)
design <- model.matrix(~clF)

library(edgeR)
d <- DGEList(counts)
d <- calcNormFactors(d)
d <- estimateDisp(d, design)
fit <- glmFit(d, design)
lrt <- glmLRT(fit, coef="clF2") #cluster 2 vs cluster 1 (end clusters)
deGenesEdgeR <- rownames(lrt$table)[p.adjust(lrt$table$PValue,"fdr")<=0.05]
mean(falseGenes%in%deGenesEdgeR) #TPR
## [1] 0.7964072
mean(deGenesEdgeR%in%nullGenes) #FDR
## [1] 0.007462687

compare GAM end point vs edgeR analysis

library(rafalib)
mypar(mfrow=c(1,1))

## compare
edgeR <- p.adjust(lrt$table$PValue,"fdr")<=0.05
tradeR <- p.adjust(endRes$pval,"fdr")<=0.05
vennC <- cbind(edgeR,tradeR)
vennDiagram(vennC)

Test for different expression pattern

resPattern <- patternTest(gamList)
hist(resPattern$pval)

deGenesPattern <- rownames(counts)[which(p.adjust(resPattern$pval,"fdr")<=0.05)]
length(deGenesPattern)
## [1] 277
mean(falseGenes%in%deGenesPattern) #TPR
## [1] 0.5528942
mean(!deGenesPattern%in%falseGenes) #FDR
## [1] 0

Monocle BEAM analysis

### old monocle BEAM analysis
library(monocle,lib.loc="~/Dropbox/PhD/Research/singleCell/trajectoryInference/trajectoryDE/methodsPaper/")
## Loading required package: Matrix
## Loading required package: ggplot2
## Loading required package: VGAM
## Loading required package: stats4
## Loading required package: splines
## 
## Attaching package: 'VGAM'
## The following object is masked from 'package:mgcv':
## 
##     s
## Loading required package: DDRTree
## Loading required package: irlba
featureInfo <- data.frame(gene_short_name=rownames(counts))
rownames(featureInfo) <- rownames(counts)
fd <- new("AnnotatedDataFrame", featureInfo)
cds <- newCellDataSet(cellData=counts, featureData=fd, expressionFamily=negbinomial.size())
cds <- estimateSizeFactors(cds)
cds <- estimateDispersions(cds)
## Removing 154 outliers
cds <- reduceDimension(cds)#, max_components = 4, method = 'ICA')
cds <- orderCells(cds)#, num_paths=2)
plot_cell_trajectory(cds, color_by = "State")

# monocle does not find a branching point on 2 components! specify 4.
BEAM_res <- BEAM(cds,  cores = 1)
## Warning in if (progenitor_method == "duplicate") {: the condition has
## length > 1 and only the first element will be used
## Warning in if (progenitor_method == "sequential_split") {: the condition
## has length > 1 and only the first element will be used
sum(BEAM_res$qval<0.05)
## [1] 404

tradeR downstream of Monocle 2

phenoData(cds)$group_id <- gid$group_id
plot_cell_trajectory(cds, color_by = "group_id")

plot(rd, col = as.numeric(as.factor(gid$group_id))+1, pch=16, asp = 1) ; legend("topleft", paste0("M",1:4), col=2:5,pch=16)

# Milestone 2 and 4 are the two branches to be compared.
ptMon <- matrix(phenoData(cds)$Pseudotime, nrow=ncol(counts), ncol=2, byrow=FALSE)
state <- phenoData(cds)$State
cellWeightsMon <- matrix(0, nrow=ncol(counts), ncol=2)
cellWeightsMon[state==1,] <- c(1/2,1/2)
cellWeightsMon[state==3,1] <- 1
cellWeightsMon[state==2,2] <- 1

gamListMon <- fitGAM(counts, pseudotime=ptMon, cellWeights=cellWeightsMon, verbose=TRUE)
## 
  |                                                                       
  |                                                                 |   0%
  |                                                                       
  |                                                                 |   1%
  |                                                                       
  |=                                                                |   1%
  |                                                                       
  |=                                                                |   2%
  |                                                                       
  |==                                                               |   2%
  |                                                                       
  |==                                                               |   3%
  |                                                                       
  |==                                                               |   4%
  |                                                                       
  |===                                                              |   4%
  |                                                                       
  |===                                                              |   5%
  |                                                                       
  |====                                                             |   5%
  |                                                                       
  |====                                                             |   6%
  |                                                                       
  |====                                                             |   7%
  |                                                                       
  |=====                                                            |   7%
  |                                                                       
  |=====                                                            |   8%
  |                                                                       
  |======                                                           |   8%
  |                                                                       
  |======                                                           |   9%
  |                                                                       
  |======                                                           |  10%
  |                                                                       
  |=======                                                          |  10%
  |                                                                       
  |=======                                                          |  11%
  |                                                                       
  |=======                                                          |  12%
  |                                                                       
  |========                                                         |  12%
  |                                                                       
  |========                                                         |  13%
  |                                                                       
  |=========                                                        |  13%
  |                                                                       
  |=========                                                        |  14%
  |                                                                       
  |=========                                                        |  15%
  |                                                                       
  |==========                                                       |  15%
  |                                                                       
  |==========                                                       |  16%
  |                                                                       
  |===========                                                      |  16%
  |                                                                       
  |===========                                                      |  17%
  |                                                                       
  |===========                                                      |  18%
  |                                                                       
  |============                                                     |  18%
  |                                                                       
  |============                                                     |  19%
  |                                                                       
  |=============                                                    |  19%
  |                                                                       
  |=============                                                    |  20%
  |                                                                       
  |=============                                                    |  21%
  |                                                                       
  |==============                                                   |  21%
  |                                                                       
  |==============                                                   |  22%
  |                                                                       
  |===============                                                  |  22%
  |                                                                       
  |===============                                                  |  23%
  |                                                                       
  |===============                                                  |  24%
  |                                                                       
  |================                                                 |  24%
  |                                                                       
  |================                                                 |  25%
  |                                                                       
  |=================                                                |  25%
  |                                                                       
  |=================                                                |  26%
  |                                                                       
  |=================                                                |  27%
  |                                                                       
  |==================                                               |  27%
  |                                                                       
  |==================                                               |  28%
  |                                                                       
  |===================                                              |  28%
  |                                                                       
  |===================                                              |  29%
  |                                                                       
  |===================                                              |  30%
  |                                                                       
  |====================                                             |  30%
  |                                                                       
  |====================                                             |  31%
  |                                                                       
  |====================                                             |  32%
  |                                                                       
  |=====================                                            |  32%
  |                                                                       
  |=====================                                            |  33%
  |                                                                       
  |======================                                           |  33%
  |                                                                       
  |======================                                           |  34%
  |                                                                       
  |======================                                           |  35%
  |                                                                       
  |=======================                                          |  35%
  |                                                                       
  |=======================                                          |  36%
  |                                                                       
  |========================                                         |  36%
  |                                                                       
  |========================                                         |  37%
  |                                                                       
  |========================                                         |  38%
  |                                                                       
  |=========================                                        |  38%
  |                                                                       
  |=========================                                        |  39%
  |                                                                       
  |==========================                                       |  39%
  |                                                                       
  |==========================                                       |  40%
  |                                                                       
  |==========================                                       |  41%
  |                                                                       
  |===========================                                      |  41%
  |                                                                       
  |===========================                                      |  42%
  |                                                                       
  |============================                                     |  42%
  |                                                                       
  |============================                                     |  43%
  |                                                                       
  |============================                                     |  44%
  |                                                                       
  |=============================                                    |  44%
  |                                                                       
  |=============================                                    |  45%
  |                                                                       
  |==============================                                   |  45%
  |                                                                       
  |==============================                                   |  46%
  |                                                                       
  |==============================                                   |  47%
  |                                                                       
  |===============================                                  |  47%
  |                                                                       
  |===============================                                  |  48%
  |                                                                       
  |================================                                 |  49%
  |                                                                       
  |================================                                 |  50%
  |                                                                       
  |=================================                                |  50%
  |                                                                       
  |=================================                                |  51%
  |                                                                       
  |==================================                               |  52%
  |                                                                       
  |==================================                               |  53%
  |                                                                       
  |===================================                              |  53%
  |                                                                       
  |===================================                              |  54%
  |                                                                       
  |===================================                              |  55%
  |                                                                       
  |====================================                             |  55%
  |                                                                       
  |====================================                             |  56%
  |                                                                       
  |=====================================                            |  56%
  |                                                                       
  |=====================================                            |  57%
  |                                                                       
  |=====================================                            |  58%
  |                                                                       
  |======================================                           |  58%
  |                                                                       
  |======================================                           |  59%
  |                                                                       
  |=======================================                          |  59%
  |                                                                       
  |=======================================                          |  60%
  |                                                                       
  |=======================================                          |  61%
  |                                                                       
  |========================================                         |  61%
  |                                                                       
  |========================================                         |  62%
  |                                                                       
  |=========================================                        |  62%
  |                                                                       
  |=========================================                        |  63%
  |                                                                       
  |=========================================                        |  64%
  |                                                                       
  |==========================================                       |  64%
  |                                                                       
  |==========================================                       |  65%
  |                                                                       
  |===========================================                      |  65%
  |                                                                       
  |===========================================                      |  66%
  |                                                                       
  |===========================================                      |  67%
  |                                                                       
  |============================================                     |  67%
  |                                                                       
  |============================================                     |  68%
  |                                                                       
  |=============================================                    |  68%
  |                                                                       
  |=============================================                    |  69%
  |                                                                       
  |=============================================                    |  70%
  |                                                                       
  |==============================================                   |  70%
  |                                                                       
  |==============================================                   |  71%
  |                                                                       
  |==============================================                   |  72%
  |                                                                       
  |===============================================                  |  72%
  |                                                                       
  |===============================================                  |  73%
  |                                                                       
  |================================================                 |  73%
  |                                                                       
  |================================================                 |  74%
  |                                                                       
  |================================================                 |  75%
  |                                                                       
  |=================================================                |  75%
  |                                                                       
  |=================================================                |  76%
  |                                                                       
  |==================================================               |  76%
  |                                                                       
  |==================================================               |  77%
  |                                                                       
  |==================================================               |  78%
  |                                                                       
  |===================================================              |  78%
  |                                                                       
  |===================================================              |  79%
  |                                                                       
  |====================================================             |  79%
  |                                                                       
  |====================================================             |  80%
  |                                                                       
  |====================================================             |  81%
  |                                                                       
  |=====================================================            |  81%
  |                                                                       
  |=====================================================            |  82%
  |                                                                       
  |======================================================           |  82%
  |                                                                       
  |======================================================           |  83%
  |                                                                       
  |======================================================           |  84%
  |                                                                       
  |=======================================================          |  84%
  |                                                                       
  |=======================================================          |  85%
  |                                                                       
  |========================================================         |  85%
  |                                                                       
  |========================================================         |  86%
  |                                                                       
  |========================================================         |  87%
  |                                                                       
  |=========================================================        |  87%
  |                                                                       
  |=========================================================        |  88%
  |                                                                       
  |==========================================================       |  88%
  |                                                                       
  |==========================================================       |  89%
  |                                                                       
  |==========================================================       |  90%
  |                                                                       
  |===========================================================      |  90%
  |                                                                       
  |===========================================================      |  91%
  |                                                                       
  |===========================================================      |  92%
  |                                                                       
  |============================================================     |  92%
  |                                                                       
  |============================================================     |  93%
  |                                                                       
  |=============================================================    |  93%
  |                                                                       
  |=============================================================    |  94%
  |                                                                       
  |=============================================================    |  95%
  |                                                                       
  |==============================================================   |  95%
  |                                                                       
  |==============================================================   |  96%
  |                                                                       
  |===============================================================  |  96%
  |                                                                       
  |===============================================================  |  97%
  |                                                                       
  |===============================================================  |  98%
  |                                                                       
  |================================================================ |  98%
  |                                                                       
  |================================================================ |  99%
  |                                                                       
  |=================================================================|  99%
  |                                                                       
  |=================================================================| 100%
resPatternMon <- patternTest(gamListMon)
resEndMon <- diffEndTest(gamListMon)

tradeR on true pseudotime

### tradeR on true pseudotime
pst <- matrix(truePseudotime, nrow=ncol(counts),ncol=2, byrow=FALSE)
gamListTrueTime <- fitGAM(counts, pseudotime=pst, cellWeights=slingCurveWeights(crv), verbose=FALSE)
# end point Test
waldEndPointResTrueTime <- diffEndTest(gamListTrueTime)
padjWaldTrueTime <- p.adjust(waldEndPointResTrueTime$pvalue,"fdr")
sum(padjWaldTrueTime<=0.05)
## [1] 297
# pattern test
patternResTrueTime <- patternTest(gamListTrueTime)
padjPatternTrueTime <- p.adjust(patternResTrueTime$pvalue,"fdr")
sum(padjPatternTrueTime<=0.05, na.rm=TRUE)
## [1] 219

GPfates

# export
logCpm <- edgeR::cpm(counts, prior.count=.125, log=TRUE)
sampleInfo <- data.frame(global_pseudotime=truePseudotime)
rownames(sampleInfo) <- colnames(counts)
write.table(logCpm, file="~/Dropbox/PhD/Research/singleCell/trajectoryInference/trajectoryDE/tradeRPaper/simulation/sim2_dyntoy_bifurcating_4/simDyntoyLogCpm.txt", row.names=TRUE, col.names=TRUE, quote=FALSE)
write.table(sampleInfo, file="~/Dropbox/PhD/Research/singleCell/trajectoryInference/trajectoryDE/tradeRPaper/simulation/sim2_dyntoy_bifurcating_4/sampleInfoSimDyntoy.txt", row.names=TRUE, col.names=TRUE, quote=FALSE)
# run GPfates python script
GPfatesWeights <- read.table("~/Dropbox/PhD/Research/singleCell/trajectoryInference/trajectoryDE/tradeRPaper/simulation/sim2_dyntoy_bifurcating_4/GPfatesWeights.txt", header=FALSE)
GPfatesBif <- read.table("~/Dropbox/PhD/Research/singleCell/trajectoryInference/trajectoryDE/tradeRPaper/simulation/sim2_dyntoy_bifurcating_4/GPfatesBifStats.txt", header=FALSE)
colnames(GPfatesBif) <- c("bif_ll", "amb_ll", "shuff_bif_ll", "shuff_amb_ll", "phi0_corr", "D", "shuff_D")

GPfates weights + tradeR

# based on true pseudotime
gamListGPfatesTrueTime <- fitGAM(counts, pseudotime=pst, cellWeights=GPfatesWeights, verbose=FALSE)
# end point test
waldEndPointResTrueTimeGPfates <- diffEndTest(gamListGPfatesTrueTime)
# pattern test
patternResTrueTimeGPfates <- patternTest(gamListGPfatesTrueTime)

plot the false genes

pdf("~/Dropbox/PhD/Research/singleCell/trajectoryInference/trajectoryDE/tradeRPaper/simulation/sim2_dyntoy_bifurcating_4/falseGenesTruePseudotime.pdf")
for(i in 1:length(falseGenes)){
   plotSmoothers(gamListTrueTime[[falseGenes[i]]])
}
dev.off()
## quartz_off_screen 
##                 2

plot the null genes (which are all permuted)

pdf("~/Dropbox/PhD/Research/singleCell/trajectoryInference/trajectoryDE/tradeRPaper/simulation/sim2_dyntoy_bifurcating_4/nullGenesTruePseudotime.pdf")
for(i in 1:length(nullGenes)){
   plotSmoothers(gamListTrueTime[[nullGenes[i]]])
}
dev.off()
## quartz_off_screen 
##                 2

plot example genes

mypar(mfrow=c(1,2))
i=1
plotSmoothers(gamListTrueTime[[falseGenes[i]]])
plotSmoothers(gamListTrueTime[[nullGenes[i]]])

Performance plots

########
# FDP-TPR
########
library(iCOBRA)
## 
## Attaching package: 'iCOBRA'
## The following objects are masked from 'package:BiocGenerics':
## 
##     score, score<-
library(scales)
truth <- as.data.frame(matrix(rep(0,nrow(counts)), dimnames=list(rownames(counts),"status")))
truth[falseGenes,"status"] <- 1
cols <- hue_pal()(10)
names(cols) <- c("BEAM", "GPfates", "edgeR", "tradeR_slingshot_end", "tradeR_slingshot_pattern", "tradeR_GPfates_end", "tradeR_GPfates_pattern", "tradeR_Monocle2_end", "tradeR_Monocle2_pattern")

### estimated pseudotime
pval <- data.frame( tradeR_slingshot_end=endRes$pval,
                    tradeR_slingshot_pattern=resPattern$pval,
                    BEAM=BEAM_res$pval,
                    edgeR=lrt$table$PValue,
                    tradeR_Monocle2_end=resEndMon$pvalue,
                    tradeR_Monocle2_pattern=resPatternMon$pvalue,
                      row.names=rownames(counts))
cobra <- COBRAData(pval=pval, truth=truth)
cobra <- calculate_adjp(cobra)
## Object up to date
cobra <- calculate_performance(cobra, binary_truth="status")
cobraplot <- prepare_data_for_plot(cobra, colorscheme=cols[match(sort(unique(cobra@roc$method)),names(cols))])
## Warning in define_colors(cobraperf = cobraperf, palette = colorscheme,
## facetted = facetted, : too few colors provided, 1 random colors will be
## added
plot_roc(cobraplot)

plot_fdrtprcurve(cobraplot, pointsize=0, yaxisrange=c(0.6,1))
## Warning: Removed 1231 rows containing missing values (geom_path).
## Warning: Removed 7 rows containing missing values (geom_point).

pEst <- plot_fdrtprcurve(cobraplot, pointsize=0, yaxisrange=c(0.6,1), xaxisrange=c(0,0.5))


### true pseudotime
pval <- data.frame( tradeR_slingshot_end=waldEndPointResTrueTime$pval,
                    tradeR_slingshot_pattern=patternResTrueTime$pval,
                    BEAM=BEAM_res$pval,
                    edgeR=lrt$table$PValue,
                    tradeR_GPfates_end=waldEndPointResTrueTimeGPfates$pval,
                    tradeR_GPfates_pattern=patternResTrueTimeGPfates$pval,
                      row.names=rownames(counts))
score <- data.frame(GPfates=GPfatesBif$D,
                      row.names=rownames(counts))
cobra <- COBRAData(pval=pval, truth=truth, score=score)
cobra <- calculate_adjp(cobra)
## Object up to date
cobra <- calculate_performance(cobra, binary_truth="status")
## column GPfates is being ignored for NBRS calculations
## column GPfates is being ignored for TPR calculations
## column GPfates is being ignored for FDR calculations
## column GPfates is being ignored for FPR calculations
cobraplot <- prepare_data_for_plot(cobra, colorscheme=cols[match(sort(unique(cobra@roc$method)),names(cols))])
## Warning in define_colors(cobraperf = cobraperf, palette = colorscheme,
## facetted = facetted, : too few colors provided, 1 random colors will be
## added
plot_roc(cobraplot)

plot_fdrtprcurve(cobraplot, pointsize=0, yaxisrange=c(0.6,1), xaxisrange=c(0,0.5))
## Warning: Removed 9811 rows containing missing values (geom_path).
## Warning: Removed 10 rows containing missing values (geom_point).

pTrue <- plot_fdrtprcurve(cobraplot, pointsize=0, yaxisrange=c(0.6,1), xaxisrange=c(0,0.5))
#
library(cowplot)
## Warning: package 'cowplot' was built under R version 3.5.2
## 
## Attaching package: 'cowplot'
## The following object is masked from 'package:ggplot2':
## 
##     ggsave
prow <- plot_grid( pEst + theme(legend.position="none") + xlab("FDP") + ggtitle("Based on estimated pseudotime"),
           pTrue + theme(legend.position="none") + xlab("FDP") + ggtitle("Based on true pseudotime"),
           align = 'vh',
           labels = c("a", "b"),
           hjust = -1,
           nrow = 1
           )
## Warning: Removed 7542 rows containing missing values (geom_path).
## Warning: Removed 7 rows containing missing values (geom_point).
## Warning: Removed 9811 rows containing missing values (geom_path).
## Warning: Removed 10 rows containing missing values (geom_point).
legend_a <- get_legend(pTrue + theme(legend.position="bottom"))
## Warning: Removed 9811 rows containing missing values (geom_path).

## Warning: Removed 10 rows containing missing values (geom_point).
p <- plot_grid( prow, legend_a, ncol = 1, rel_heights = c(1, .2))
p

# png("~/Dropbox/PhD/Research/singleCell/trajectoryInference/trajectoryDE/tradeRPaper/simulation/sim2sim2_dyntoy_bifurcating_4/performance_sim2Dyntoy.png", width=7,height=6, units="in", res=300)
# p
# dev.off()

#composite one plot
pval <- data.frame( tradeR_slingshot_end=endRes$pval,
                    tradeR_slingshot_pattern=resPattern$pval,
                    BEAM=BEAM_res$pval,
                    edgeR=lrt$table$PValue,
                    tradeR_GPfates_end=waldEndPointResTrueTimeGPfates$pval,
                    tradeR_GPfates_pattern=patternResTrueTimeGPfates$pval,
                    tradeR_Monocle2_end=resEndMon$pvalue,
                    tradeR_Monocle2_pattern=resPatternMon$pvalue,
                      row.names=rownames(counts))
score <- data.frame(GPfates=GPfatesBif$D,
                      row.names=rownames(counts))
cobra <- COBRAData(pval=pval, truth=truth, score=score)
saveRDS(cobra,"~/Dropbox/PhD/Research/singleCell/trajectoryInference/trajectoryDE/tradeRPaper/simulation/sim2_dyntoy_bifurcating_4/cobra.rds")
cobra <- calculate_adjp(cobra)
## Object up to date
 cobra <- calculate_performance(cobra, binary_truth="status")
## column GPfates is being ignored for NBRS calculations
## column GPfates is being ignored for TPR calculations
## column GPfates is being ignored for FDR calculations
## column GPfates is being ignored for FPR calculations
 cobraplot <- prepare_data_for_plot(cobra, colorscheme=cols[match(sort(unique(cobra@roc$method)),names(cols))])
## Warning in define_colors(cobraperf = cobraperf, palette = colorscheme,
## facetted = facetted, : too few colors provided, 1 random colors will be
## added
  plot_roc(cobraplot)

 plot_fdrtprcurve(cobraplot, pointsize=0, yaxisrange=c(0.8,1), xaxisrange=c(0,0.5))
## Warning: Removed 13269 rows containing missing values (geom_path).
## Warning: Removed 22 rows containing missing values (geom_point).

# pAll <- plot_fdrtprcurve(cobraplot, pointsize=3/2, yaxisrange=c(0.8,1), xaxisrange=c(0,0.5))
# png("~/Dropbox/PhD/Research/singleCell/trajectoryInference/trajectoryDE/plots/performanceSim1All.png", width=7,height=6, units="in", res=300)
# pAll
# dev.off()